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ABSTRACT 



The practical and rigorous solutions of the potential 
problem associated with the harmonic oscillation of a sub- 
merged rigid body .of arbitrary shape is presented. The use 
of Green's function reduces the determination of the poten- 
tial to the solution of an integral equation. The integral 
equation is solved numerically and the dependency of the 
hydrodynamic quantities such as added mass and damping 
coefficients of the object on the frequency of the oscilla- 
tion is established. 

Several checks are made on the numerical results. These 
include the Haskind's relations check, an energy check for 
the radiation problem, and comparisons with closed form 
solutions. All the checks and comparisons are successful and 
it appears that the numerical procedure employed yields valid 
and accurate results. 
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NOMENCLATURE 



Symbol 


Description 


a 


characteristic length of object 


a 


wavelength parameter, a = 2ira/L, dimensionless 


C. 

1 


th 

wave force or moment coefficient in the l 
direction, dimensionless 


di 


elementary dimensional surface area 


ds 


elementary dimensionless surface area 


E. 

1 


energy transmitted over one period 


v f 

F. 

l 


distribution function, dimensionless 

t b 

wave force or moment component in the l 
direction, dimensionless 


g 


acceleration of gravity 


G 


Green’s function 


h 


water depth 


h 


relative depth of water, dimensionless 


h. 

3 

H 


functions defined by (19) , dimensionless 
wave height 


Im 


imaginary part 


n x ,n y ,n z 
P' 


unit normal vectors, dimensionless 

dynamic pressure associated with the wave 
interaction with the fixed object 


P. 

3 


"t h 

dynamic pressure associated with the i 
mode of oscillation, dimensionless 


q 


fluid velocity vector 


Re 


real part 


r 


plan polar coordinate 


r 


plane polar coordinate, dimensionless 
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Symbol 

r. 



s 

s 

t 

T 



u . 

1 

X. 0 
1 

X.° 

l 

x,y,z 

x,y,z 







N. . 
11 

6 



€,t],C 

0 



9 . 
1 



V 



p l 



Description 

horizontal distance between points (x,y,z) 
and dimensionless 

surface area of the object 

surface area of the object, dimensionless 

time 

wave period 

velocity potential for incident wave, dimensionless 

velocity potentials for the radiation and scatter 
problems, dimensionless 

amplitude of motion in i L direction 

amplitude of linear or angular motion, dimensionless 

spacial variables 

spacial variables, dimensionless 
t ll 

i r component of force or moment coefficient 
associated with the jth component of oscillation 
of the object 

t h 

i component of the dynamic force or moment 
due to the j th mode of oscillation 

added mass coefficient 

added mass coefficiemt, dimensionless 

damping coefficient 

damping coefficient, dimensionless 

Dirac delta function 

rectangular Cartesian coordinates, dimensionless 
plane polar coordinate 

angular displacement of oscillating object 
a a/g = a tanh (ah) 



2 2 1/2 

cr + v) 
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Symbol 



Description 

a angular frequency, a = 2tt/T 

<f> o velocity potential for incident wave 

<J> . velocity potential for oscillation of object 

- 1 or scatter 

total velocity potential for diffraction problem 



Subscripts 



i,j used to denote either the mode or oscillation 

Cor scatter) or the direction of a force or 
moment component, depending on context 
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I. INTRODUCTION 



The need for the utilization of the natural resources 
imbedded in coastal waters has led to the need for greater 
knowledge of the interaction of surface gravity waves with 
large submerged objects such as oil storage tanks and pipe- 
lines. A new problem arises in the case of large bodies 
immersed in shallow water, the scattering of the incident 
wave and the free surface effect become an important considera- 
tion in the calculation of wave forces on fixed objects and 
hydrodynamic forces produced by oscillation of the body in 
otherwise still water. 

The analysis which considers these two factors has come 
to be known as diffraction theory. In this approach separa- 
tion and viscous effects may be neglected, as was confirmed 
by Sarpkaya and Garrison [1], provided the amplitude of the 
motion is sufficiently small. The problem is set up in 
terms of a velocity potential and the potential function is 
determined by use of a source distribution. The pressure 
on the immersed surface can be found once the source strengths 
are known by use of Bernoulli’s equation and the resulting 
forces and moments are determined from the pressure distribu- 
tion by surface integration. 

Recently, Garrison and Seetharama Rao [2] worked out an 
analysis for wave interaction with submerged objects based on 
linear diffraction theory and applied the analysis to calculate 
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wave forces acting on a submerged hemisphere. They showed 
that for small values of the wave height to sphere diameter 
ratio, the wave force coefficients can be well represented 
by the diffraction theory and correlated as functions of the 
parameters, 2TTa/L and h/a where h denotes the water depth, a 
the sphere radius, and L the wave length. 

Aside from the wave force problem, an additional problem 
may be considered which is referred to as the ’’radiation 
problem." This involves the effect of a rigid body in calm 
water undergoing forced harmonic oscillations in its various 
degrees of freedom. Of interest here are the hydrodynamic 
forces which may be characterized by the coefficients of 
added mass, added moment of inertia, and damping. Although 
often solved independently, the diffraction and radiation 
problems may be solved concurrently because they are mathe- 
matically similar. 

Relations which eliminate the need to solve the problem 
of wave diffraction past an object by solving the forced- 
oscillation problem in still water were developed by 
M. D. Haskind in 1957. They relate the force and moment 
components of the diffraction problem to the asymptotic 
velocity potentials for the corresponding radiation problem. 
That is, if the waves produced by the oscillation of an 
object in still water are known, Haskind's relations may be 
used to obtain the wave forces (or moments) produced by 
surface waves interacting with the same object held fixed. 
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Conservation of energy dictates that a balance must exist 
between the energy required to oscillate an object, and the 
wave energy transmitted by surface waves across some control 
volume surrounding the object but at a large radial distance. 

The energy required to oscillate the object is directly 
related to the damping coefficient, and the damping coefficient 
may be obtained by integration of the pressure over the immersed 
surface. In view of conservation of energy this same energy 
appears at a large distance from the object in the form of 
radiated wave energy. Therefore, the damping coefficient may 
be obtained by accounting for this energy flux. The compari- 
son of this coefficient obtained from these two methods, i.e., 
the near field and far field solutions, provides a valuable 
check on the accuracy of the numerical results generated. 

The application of Haskind's relations and energy considera- 
tions to calculate wave forces on large submerged objects is 
the purpose of this thesis. A computer program was developed 
to carry out the numerical calculations and to compare the 
results with those from diffraction and radiation theory. 
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II. THEORETICAL CONSIDERATIONS 



A. A REVIEW OF WAVE FORCES ON SUBMERGED OBJECTS 

Consider a rigid body of arbitrary shape having a charac- 
teristic size a submerged in water of depth h as illustrated 
in Fig. 1. It is assumed that the amplitude of the motion 
of the fluid is sufficiently small in relation to the size, 
a, that viscous effects can be neglected. The object may 
or may not be in contact with the free surface or the bottom. 
The problem set forth herein deals with the fluid motion and 
forces induced by the small amplitude oscillation of the 
object in its six degrees of freedom as well as the fluid 
motion associated with the interaction of the fixed object 
with a train of regular waves. 

The small amplitude oscillatory motion of the body with 
frequency about its equilibrium position is described by the 
relationships 

X ± (t ) = X° Re [e" 1CTt ] , i = 1, 2, 3 (la) 

9 i Ct ) = 9° Re[e‘ iat ] , i = 4, 5, 6 (lb) 

__ O 

where denotes the amplitude of motion in translational 

O 

oscillation and 0^ denotes the amplitude of motion in rota- 
tional oscillation. The subscript i = 1, 2, 3 corresponds to 
oscillation in the x, y, z direction, respectively, and 
i = 4, 5, 6 corresponds to angular oscillation about axes 
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X 




15 



Figure 1: Definition Sketch 



passing through the origin of the object parallel to the 
x, y, z axes, respectively. The second problem considered 
here is the interaction of a train of regular surface waves 
with the object fixed in space. The incident waves of wave 
height H and wave length L are assumed to progress in the 
positive x-direction and the fluid motion is assumed to be 
incompressible, irrotational , and harmonic with time dependence 
e 10t in all cases. Consequently, a velocity potential exists 
such that the fluid velocity vector may be defined as 

qj = Re[V<j) i e' iot ] , i = 1, 2, ...6 (2) 

where <j>^ denotes the velocity potential associated with the 
motion induced by oscillations in the six degrees of freedom. 

In the case of regular wave interaction with the fixed object, 
the velocity potential is given by the sum 

<P’ = <P 0 + <t>y (3) 

where <f> denotes the velocity potential of the incident wave 
in the absence of the object and <f>y denotes the velocity 
potential of the scattered wave due to the presence of the 
rigid body. For this case, the fluid velocity vector is 
given by 

q' = Re[VO 0 + <J> 7 )e" iat ] . (4) 

Having assumed irrotationality of the fluid, <J>^ must satisfy 
the continuity equation which takes the form of the Laplace 
equation, V 2 <jK = o 



( 5 ) 



within the fluid region. 



Assuming that the amplitude of the motion is small, the 
velocity squared terms in Bernoulli's equations may be neg- 
lected. Therefore, the dynamic pressure is given as 

P. = Re[pa<f> •e'" icrt 3 , j =1, 2,... 6 (6) 

For the second problem, involving wave interaction with the 
fixed object, the pressure is given by 

P’ = Re[pa(<f> o + <f> 7 )e‘ icr,: ] . (7) 

The velocity potential must satisfy certain boundary 
conditions in addition to Eq. (5). These include the 
linearized free surface boundary condition, 

_ _ rr 2 

— 3. (x, 0, z) p. (x, 0, z) = 0 , j = 1, 2,. ..7, (8) 

3y S3 

as well as the kinematic boundary condition on the bottom, 

dp. 

— i (x, -h, z) = 0 , i = 1, 2, ...7. (9) 

3y 

Since 4> Q satisfies both Eqs . (8) and (9), it follows, there- 

fore, that pj must also satisfy these conditions. 

In addition, 4 k must satisfy the boundary condition on 
the surface of the body defined by s(x, y, zj =0. This con- 
dition states that for an impermeable boundary, the normal 
velocity of the fluid must just equal the velocity of the 
surface normal to itself. In the case" of oscillation of the 



rigid body, the kinematic boundary conditions on the surface 



of the object take the form 



3<f>, i o 

(x, y, z, t) = X, n 



3n 



1 x 



3 4 > 2 : o 

-T- U, y, Z, t) = X ? n 
3n y 



3 4> - 



3n 



9<f> 



3n 



(x, y, z, t) = X,n 



3 ii 

z 



(x, y, z, t) = 0 4 [(h + y)n z - xn y ] 



( 10 ) 



d<P, 



3n 



(x, y, z, t) = 9 5 [zn x - xnj 



^^6 * - - 

-T“ (x, y, z, t) = 9 fi [xn - (h + y)n ] 

5n o y x 



where n = in + in + kn denotes the unit normal vector on 
x J y z 

the surface of the object directed outward into the fluid. For 
the case of the interaction of regular waves with a fixed, 
rigid body the normal velocity must be zero which may be 
written 



3<J>____ 3<f> 

— Z~ (x> y> z, t) = - —2- (x, y, z, t) . (11) 

3n 3n 

Finally, the velocity potential must satisfy the radiation 
condition which permits only an outgoing progressive wave at 
a large distance. Hence it is required that <J)j have the 
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asymptotic form which satisfies 



- 1/2 - - 

. r ~ * i a cosh[k(h+v)l i(kr.-at) ^ n 

<p.(r, 0, y) - A-CGJr, — - — <- J - L & y 1 ' 0 as r, 

3 3 1 cosh (kh) 1 



3 = 1, 2, 



- 2 - 21/2 -- 
7 where r^ = (x + z ) ' and 0 = arctan (z/x) 



where (0) is some unknown complex function of 0. 

The wave number is defined as k = 2 tt/L where L is the wave 
length and is related to the frequency of disturbance according 
to 



^— = k tanh (kh) 



(13) 



The velocity potential of the incident wave alone progres- 
sing in the positive x-direction is given by 



4> 0 (x, y, z) = 



ign cosh[k(h + y) ] e ~ikx 



cosh (kh) 



(14) 



where q = H/2 represents the amplitude of the incident wave, 

H being the wave height. 

In carrying out the solution and to show more clearly the 
dependence of the solution on the relative wave length parame- 
ter, a = 2Tra/L = ka, and relative depth, h = h/a, a being a 
typical length of the body, it is convenient to make the space 
variables and amplitudes dimensionless, i.e., 



x = x/a , y = y/a , z = z/a , r = r/a 

X. = X./a , i = 1, 2, 3 and X. =0. , i=4, 5, 6 

l i * * ’ i l ’ ’ ’ 



( 12 ) 
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and then to introduce 

O 

ic^jU, y, z) / gaXj 

o 

io<}) o (x, y, z)/ gar) 

O 

io<t> 7 (x, y, z) / gan 



the pressure function u^ by 
= a tanh (ah) u^ (x, y, z), 
= -au Q (x, y, z) 

= -au y (x, y, z) . 



= 1 , 2 



The complex dimensionless dynamic pressure amplitude can now 
be written by use of the linearized form of Bernoulli's 
equation as 



P’ 



a tanh (ah) Uj 

cosh[a(h + y) ] 
cosh (ah) 



C x, y, z) , j = 1, 
e iaX -au ? (x, y, z) 



2 , 



.6 



where the complex amplitude of the pressure p is defined as 



-44 = Re[p • e" iot ] 

PgaXj ] 



P' 

Pgan 



Re[p'e 



The boundary value problem which describes the fluid 
motion arising from the oscillation of the rigid body in its 
six degrees of freedom as well as the scattering of the 
incident wave can now be written concisely in terms of dimen 
sionless parameters. The potential u^(x, y, z), j = 1, 2, 
continuous in the fluid region is sought such that 

V 2 Uj (x, y, z) = 0 in y < 0 (a 



« . . 6 

(15) 



(16) 



(17) 



) (1 
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3u. 

3v* 



(x, 0, Z) - a tanh (ah) (x, 0, z) = 0 (b) 



3u . 



3y 



i (x, -h, z) = 0 outside s(x, y, z) = 0 



(c) 



3u . 

1 

3n 



(x, y, z) = h. (x, y, z) on s(x, y, z) = 0 (d) 



Uj (r, 6 , y) - Dj (6) 



1/2 cosh[a(h + v) ] iar 
cosh (ah) e 



i -> 0 as r 00 (e) 



where (0) is some dimensionless complex function of 0, 
s(x, y, z) = 0 represents the surface of the object in 
its undisturbed position and iu denotes the prescribed 
function which depends on the mode of oscillation (j = 1, 2, . 

the subscript 7 corresponds to scattering of the regular inci- 
dent wave due to the presence of the rigid body. For the six 
degrees of freedom and scattering, the functions h j are given, 
respectively, by 



6 ) 



h. = n 
1 x 



, h _ = n 
’ 2 y 



h „ = n 



h 4 = (h+y)n z - zn y , h g = zn^ - xn z , h & = xn y - (h+y)n^ (19) 



y 



h 7 = c o s h T a K y t n y sinh[a(h+y)] + in x cosh[a(h+y)] e 



x 
iax 



B. REPRESENTATION OF THE POTENTIAL 



The solution to the boundary value problem (18) may be 
attained by use of a Green's function having the physical 
interpretation of a point wave source of unit strength. These 
sources are distributed over the surface of the object 
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according to the source strength function f so that the 
potential at some point x, y, z within the fluid is given by 
the surface integral 



u . 
1 




f(?» n, C) G(x, y, z; n , ?)ds 



C20) 



s 

where (£, n, C) represents coordinates on the surface of the 

object, f(£, n, £) represents a distribution or weighting 

- - 2 

function and ds = ds/a denotes the dimensionless surface area 
element. 

The Green's function is defined, therefore, as the function 
which satisfies 



V 2 G (x, y, z; £, n, C) = <5 (x - £) 5(y - n) 6(z - C) (21) 

where 6 is the Dirac delta function, as well as the boundary 
conditions (18b, c, e) . Such a function is given by Wehausen 
and Laitone [3] as 

G (x, y, z; £, n, c) = + G*(x, y, z; C, n, C) (22) 

where 



G*(x, y, z; ?, n, ?) = 



R 1 



x on V f (y + V)e _yh cosh [y (h + n)] cosh [y(h + y) ] rnridu 
2P ‘ V * / y sinh (ph) cosh (yh) J Q (yr)du 

O 



i 2?r(a 2 - v 2 ) cosh [a(h + n)] cosh [a(h + y) 1 j 

a 2 h - v 2 h + v ° 



(23) 
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R = [(x - C) 2 + (y - n) 2 + (z - c) 2 ] 1/2 

R’ = [(x - ?) 2 + (y + 2h + n) 2 + (z - c) 2 ] 1/2 

r = [(x - £) 2 + (z - c) 2 ] 1/2 

v = 2—2. s a tanh (ah) 
g 

and P.V. indicates the Cauchy principle value of the integral. 
An alternate series form of Green's function is also given 
by Wehausen and Laitone as 

~ f 2 2-. 

G(x,y,z;£,n,C) = -4^ — 2 — - cosh[a(h+y)] cosh [a (h+n) ] • 

a n-v h+v 

[Y Q (ar) - iJ Q (ar) ] (24) 

(y 2 -v 2 ) 

+ 4 I — 2 cos[y k (h+y)]cos[p k (h+n)]K (y k r) 

k=l P k h-v h+v 

where J q and Y q denote, respectively, Bessel functions of the 

first and second kind of order zero and K denotes the modi- 

o 

fied Bessel function of the second kind of order zero. The 
quantities y k are the real positive roots of the equation 

y k tan (y k h) + v = 0 . (25) 



The solution to the boundary value problem now rests on 
the determination of the source strength function f. Taking 
the normal derivative of the potential u^ and applying boundary 
condition (18d) yields the following integral equation from 
which f is to be determined: 
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= 1, 2, . . .7 (26) 



1 _ 

4tt 




£ j ^’ n ’^f^ X ’ y ’ Z; ^ ,n ’ C ) ds 



h - (x,y,z) 



j 



where 3G/3n may be obtained by straight forward differentiation 
of (22) and (24). 

C. NUMERICAL SOLUTION 



A numerical procedure can be devised by approximating the 
actual contour by a profile of a finite number of facets. 
Proceeding in this direction, the surface of the object is 
partitioned into N area elements of size As ^ where the sub- 
script takes on values j = 1, 2, ...N. Since the source 

strength function f, occurring in the integrals in (20) and 
(25), is a continuous, well-behaved function, these integrals 
may be approximated by the following summations: 



f . 




h . 



ni 



where 



a 



ij 




u . = f . 3 - • 

ni ni ij 




(27) 
C 2 8 ) 

(29) 

(30) 



In these expressions, u ^ denotes the potential at the i 
nodal point on the object associated with either the n^.^ mode 
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* 










of o. illntic.'*, or scattering (n = 7) of the incident wave. 

The Titegr -“t-ons are carried out over the finite surface 

aTea elem As. and f. denotes the value of f at the 

.. 3 J n 

central point? of the area element. The approximation made 
heTc is conciitent in that as N approaches infinity, the 
approximations given in (27) and (28) approach the exact 
forms . 

■ 'nonce th.l : r: at r: : - ces a, 8, and h are numerically evaluated, 

(27;"can c;u. iriy be inverted by use of a standard computer 

subroutine to obtain the solution for f. u is then obtained 

n 

from (29). There are, however, certain difficulties connected 
with the evaluation of (29) and (30) which must first be con- 
sidered. Specifically for the special case when i = j , the 
1/R"and 3(l/R)/3n terms are singular as R approaches 0, and 
•t sectuiuly in the evaluation of (29) and (30) due to the singu- 
lar nature ol (l/(y tanh (yh) - a tanh (ah)] at y = a occurring 
in the infinite integral in (23) and its normal derivative. 

These difficulties are, however, circumvented by subtracting 
out '.the singularity and carrying out its integration analytical- 
ly. These details are carried out in some, yet unpublished, 
notes of X... J’. Garrison. 

D. HYDRODYNAMIC FORCES AND MOMENTS 

The forces and moments caused by the dynamic fluid pressure 
acting upon the immersed surface of the object may be obtained 

f r cm tli e ■ in t ngr a 1 s , 
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...6 (31) 



F i j (t) = - (1 or a) JJ P.h.di , i, j = 1, 2, 

F ± (t) = - (1 or a) JJ P , h i di , i = 1, 2, ...6 (32) 

til 

where denotes the i component of wave force (or moments) 

til 

and F^. denotes the i L component of forces arising from the 
j*'* 1 component of body motion. The coefficient 1.0 is used for 
the case of a force (i = 1, 2, 3) whereas a is used when F 
denotes a moment (i = 4, 5, 6). 

For purposes of presentation of the numerical results, 
dimensionless force coefficients are defined as 



and 



C . 

l 



F. • , 

imax io. 

_ 3 o e 1 

Pga n 



i = 1, 2, 3 




F. 

imax 

- 4 ° 
PgaXj 






i = 1, 2, 3 

j 1, 2 , ...6. 



(33) 



(34) 



The corresponding expressions for the moment coefficient are 



C. 

l 



F. . . 

imax i<5. 

- - 4-0 e 1 
Pga n 



i = 4, 5, 6 



(35) 



and 




F. 

imax 

- 4 Y o 
Pga Xj 






) 



i = 4, 5, 6 

J *“ l y 2 , •••6* 



(36) 



til 

The complex coefficients relate to the i component of 
wave force (or moment) coefficient while denotes the i 

component of force (or moment) coefficient associated with 

t^ 

the j component of oscillation of the object. The phase 



th 
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shift angles 5^ and 5^ relate the phase of the force to the 
crest of the incident wave and displacement of the body, 
respectively. The real and imaginary parts of the dimension- 
less force coefficients, and N^. , are called the added 

mass and the damping coefficients, respectively. 

Using (31) and (32) in conjunction with (16) and (17) 
and the definitions (33-36), the force coefficients may be 
written in dimensionless form as 



C i = JJ [au 7 (x,y,z) - C °cosh '' ( ~ ih) 1 h i Cx,y,z)ds, i= 1,2, ...6 (37) 



and 



C— = ff u j (x,y ,z)h i (x,y,z)ds , i, j = 1, 2, ...6. 



(38) 



Once Uj(x,y,z) is obtained from (28), the coefficients can be 
obtained from (37) and (38) by numerical quadrature. 

E. HASKIND 1 S RELATIONS AND ENERGY CHECK 

Although it may be presupposed that the numerical solution 
outlined here will converge upon increasing the number of 
partitions, modern computers are still limited by storage 
capacity and computer time is expensive. In view of this, 
it is necessary to keep the partition size as large, i.e., 
the number of partitions as small, as accuracy considerations 

permit. The determination of the effect of the partition 
size on accuracy becomes, therefore, an important considera- 
tion in order that practical limits may be established. One 
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available method is the comparison of the numerical results 
with analytical results where closed form solutions exist. 
This approach, however, is limited to a few simple shapes; 
for more general configurations no such check, of course, 
exists . 

A method which is valid for checking the soundness of 
the numerical results of any shape is the use of the so- 
called Haskind's relations as well as an energy balance, both 
of which are developed in the Appendix. The former is used 
to verify the solution for the diffraction problem and 
the latter the solution of the radiation problem. Using the 
asymptotic form of Green's function as given in (24) in con- 
junction with (20) the following relationship for the damping 
coefficient is so obtained 



the damping coefficient to the far field behavior of the solu- 
tion. It is in effect, an energy balance between the energy 
required to oscillate the body and the energy transmitted 
across some control surface a large distance from the body. 

A relationship somewhat similar to (39), known as Haskind's 
relations, may be obtained for the wave force (or moment) 




cosh [ a (h+q) ] e 



-ia Pl cos(3-9j ds 2 d0 



(39) 



2 7 1/2 

where = (£ + £ ) ' and 3 = arctan (£/£). This relates 
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coefficient and relates the wave produced at infinity by the 

th t h 

body oscillating in the j mode to the j component of 
wave force as 



C j cosh (ah) JJ f.(€,n,C) 

s J 

cosh [a(h+n) ]e" iap l cos ^ -Tr -*ds. 



(40) 



This result represents the Haskind’s relation for finite depth 

as given by Seetharama Roa [4] . 

Equations (39) and (40) represent relations for the 

damping and wave force (or moment) coefficient based on the 

behavior of the far field solution. A comparison of these 

results with N. . and C. obtained from an integration of the 
ij 3 

pressure over the immersed surface, i.e., as obtained from 
(37) and (38) , provides a convenient and valuable self-check 
on the accuracy of the numerical results. These results are 
not limited to special configurations and may be applied to 
arbitrary shapes. Equation (40) is, however, limited to 
symmetry with respect to the x-y plane, a condition which is 
generally satisfied by practical shapes. 



III. COMPUTER PROGRAM 



In order to feasibly utilize the method described in 
section D of THEORETICAL CONSIDERATIONS , a computer program 
was developed to carry out the indicated calculations. Although 
accuracy was a controlling factor, as it should be, the pri- 
mary requirement was to hold the computer run time to a 
reasonable value. To achieve this, every advantage was taken 
to employ the symmetry in the equations and to eliminate redun- 
dant computations by generating certain constants and matrices. 

The program consists of several subroutines, each solving 
a specific part of the problem. This allowed the development 
of the program to be simplified and the use of it was made 
easier in that only a small main program was required to 
execute the program. Subroutine GEODAT reads the input geo- 
metrical data, generates the matrix h as defined in (19), and, 
in conjunction with GCOEFF, calculates certain geometrical 
parameters. GARRIS, the primary subroutine, then generates 
the a and 3 matrices by calling either of two subroutines. 

The first (GREEN) calculates G and 3G/9n based on the integral 
form given in (22) and (23) whereas the second (GREENS) makes 
evaluations on the basis of the alternate series formulations 
given in (24). For elements of the a and 8 matrices corres- 
ponding to small values of (ar) GREEN is utilized, while for 
large values of (ar) GREENS is used. With the exception of 
the diagonal elements, the majority of the elements of a and 




Figure 2: Program Flow Chart 
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8 are calculated by the latter subroutine. This situation 
proved to be most fortunate because the series form converges 
rather rapidly and, consequently, requires much less computer 
time than the integral form. 

In order to check the validity of the two subroutines, 
values of G were computed. A comparison between GREEN and 
GREENS, Figs. 3, 4, S f for constant a and decreasing r, shows 
that at small values of a, the two methods give similar 
results whereas at intermediate values of a, the GREEN solu- 
tion becomes very erratic at larger values of r. At high 
values of a, however, the two subroutines render similar 
solutions despite the higher frequency of oscillation. This 
stems from the fact that the amplitude of oscillation decays 
extremely fast. 

Subroutine COMAT then determines the source strength 
function f by solving the complex matrix equation (27). 

Once the matrix f has been generated, the potential function 
u is calculated by (28). The pressure at points on the 
immersed surface as well as the resulting forces and moments 
are then determined in subroutine FORCES. 

Finally, the Haskind’s relations and the energy check are 
employed to calculate the wave force and moment coefficients 
and the damping coefficients in COCHCK. Here also, a compari- 
son is made between these values and those found in FORCES. 
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Figure 3: Comparison of Equations (22) and (24); a = 0.2 
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Figure 4: Comparison of Equations (22) and (24); a = 1.0 
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Figure 5 : Comparison o£ Equations ( 22 ) and (24) for Large a. 
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IV. DISCUSSION OF RESULTS 



A. VERTICAL CYLINDER 

In order to study the accuracy of the computer program, 
numerical results were generated for both the wave force 
problem and the problem involved in the oscillation of the 
object in otherwise still water. Simple geometric shapes, 
including vertical circular cylinders and submerged spheres, 
were studied and compared with closed form solutions where 
they exist. 

As indicated in Section II, the accuracy of the numeri- 
cal results is dependent on the size of the grid. Garrison 
[Ref. 5], for one, found that the fineness of the grid size, 
N, must be increased as a increases. This results from the 
fact that the kernels of the integral equations oscillate 
rapidly as the parameter a increases. 

A series of computer runs was made for a vertical circu- 
lar cylinder resting on the bottom and piercing the free 
surface. Two runs consisting of cylinders with grid sizes 
N = 120 and N = 252 were made for a relative depth of h = 2.0 
in order to determine the effect of the grid size. Another 
was made with N = 120 and h = 6.0 to determine the effect 
of depth. A fourth run was made with N = 240 and h = 4.0, 
but this time the grid was divided into f- c 'izes. The 
bottom half of the cylinder had N = 96 with an element 
aspect ratio (width/height) of 0.525 whereas the top half had 
N = 144 and an aspect ratio of 0.785.. 
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Figure 6: Definition Sketch 
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Within the range of practicality of a, a = 0 to a = 3.0, 
the numerical results for the horizontal wave force coeffi- 
cient and the coefficient of moment about a bottom line axis 
compared favorably with the closed form solution given by 
MacCamy and Fuchs [Ref 6], Figs. 7 and 8. These results 
show that within this range of the parameter a, there is such 
a small deviation in the values obtained for the different 
grid sizes that the two curves are superimposed. As the 
frequency is increased, however, it was found that deviation 
from the closed form results increased more rapidly for 
the coarser grid than the finer one. It is also note- 
worthy that the magnitude of the moment coefficient relative 
to the wave force coefficient is small at low frequencies 
and increases as the frequency increases. This result is 
as expected since at small values of wavelength, wave action 
is concentrated near the free surface whereas at large 
values of wavelength, wave action is felt at deeper depths, 
thus decreasing the moment arm. 

Figures 9 and 10 represent the hydrodynamic force and 
moment coefficients associated with the cylinder oscillating 
in still water. These show that for zero frequency, the 
force coefficient, C^, approaches approximately 2ir and 4ir 
for the h = 2.0 and h = 4.0 cylinder, respectively. At 
h = 6.0 a value of approximately 6tt was obtained. When this 
is converted to an added mass coefficient the result is 1.0. 
This is as expected since at zero frequency the free surface 
acts as a rigid boundary and the flow corresponds to 
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Figure 7: HORIZONTAL WAVE FORCE COEFFICIENT FOR A VERTICAL CIRCULAR CYLINDER 
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Figure 8: HORIZONTAL WAVE MOMENT COEFFICIENT FOR A VERTICAL CIRCULAR CYLINDER 
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Figure 9: HYDRODYNAMIC FORCE COEFFICIENT FOR 

A VERTICAL CIRCULAR CYLINDER 
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HYDRODYNAMIC MOMENT COEFFICIENT 




Figure 10: HYDRODYNAMIC MOMENT COEFFICIENT FOR A VERTICAL 

CIRCULAR CYLINDER 
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two-dimensional flow past a circular cylinder. The classical 
value for the added mass coefficient in this case is 1.0. 

B. SPHEROID 

A series of runs for varying a was made for a half 
submerged sphere. One consisted of a grid size of 200 
points with the relative depth h = 5.0 and the others were 
made with a grid size of 264 points and h = 5.0 and 
h = 10.0. Havelock [Ref. 7] determined the added mass and 
damping coefficients for a floating, heaving sphere. His 
results were later confirmed by Kim [Ref. 8] who also 
obtained such values for floating, surging (or swaying) 
spheroids. Although the results found by Havelock and Kim 
were for infinite depth, the values found in this study for 
h = 5.0 and 10.0, Figs. 11-16, compare quite favorably. 

There is apparently little depth effect beyond h = 5.0 
except at small values of a where it is noted the curves 
corresponding to h = 5.0 and h = 10.0 differ. 

Kim [Ref. 9] also determined the horizontal and vertical 
wave force coefficients as did Milgram and Halkyard [Ref. 10]. 
Comparison of these values with those calculated in this 
study, Figs. 15 and 16, show good agreement with the former, 
but some difference was noted with the latter, especially 
for the vertical force coefficient where the value of 
Milgram and Halkyard was almost three times greater at a = 2.0. 
It would appear the Milgram and Halkyard results are in error. 
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Figure 11: COF FFICIF.NT OF DAMPING FOR A IIF.AVING SPIIFRF. 
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Figure 12: COEFFICIENT OF ADDED MASS FOR A HEAVING SPHERE 
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Figure 13: COEFFICIENT OF DAMPING FOR A SURGING SPHERE 
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Figure 14: ADDED MASS COEFFICIENT IN SURGE FOR A FLOATING SPHERE 
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Figure 15: HORIZONTAL WAVE FORCE COEFFICIENT FOR A FLOATING SPHERE 
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Figure 16: VERTICAL WAVE FORCE COEFFICIENT FOR A FLOATING SPHERE 



For a fully submerged spheroid there were only two 
available checks known. One is that the value of the added 
mass coefficient for a sphere in an infinite fluid is 0.5. 
Figure 18 shows that for the numerical procedure used in 
this study a value of 0.51 was obtained for a relative depth 
of 5.0. The second test involves the use of Haskind's rela- 
tions which relates the vertical force on an object fixed 
in a monochromatic wave train with the heave damping 
coefficient such that 



N * 2 



a 

2 



s inh (2ah'> 



2ah + sinh (2alTJ 2 



(41) 



where denotes the damping coefficient in heave. A 

comparison of the damping coefficient found in this study 
with that using (41) is shown in Table 1. Using the results 
of the diffraction theory as the standard, the percentage 
deviations from it for the Haskind's relations are also 
shown. The comparison shows that for values of a greater 
than 0.5, all of the results obtained were within 4% of the 
taken standard. For values of a less than 0.5, the differ- 
ences were considerable and this was believed to be due to 
the low value of the relative depth h. A similar compari- 
son was made for the 264 point floating sphere using h = 10.0. 
The results are shown in Table 2. Since the present results 
satisfy Haskind's relations, t may be concluded the Milgra. 
and Halkyard's results are definitely in error. 
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a TANH (ah) 

Figure 18: HYDRODYNAMIC FORCE COEFFICIENTS FOR A SUBMERGED SPHERE; 200 PTS 



Table 1: Damping Coefficient in Heave for a Submerged Sphere 



N = 200 h = 5.0 



a 


N 2 2 


N* 

2 2 


% deviation 


0.20 


0.0267 


0.0199 


25.46 


0.30 


0.1033 


0.0914 


11.52 


0.50 


0.3382 


0.3252 


3.84 


0.80 


0.4433 


0.4329 


2.34 


1.00 


0.3745 


0.3698 


1.25 


1.50 


0.1610 


0.1589 


1.30 


2.00 


0.0514 


0.0515 


-0.19 


3.00 


0.0033 


0.0034 


-3.03 


4.00 


0.0002 


0.0002 


0.0 


5.00 


0.0 


0.0 


0.0 



Table 2: Damping Coefficient in Heave for a Floating Sphere 

N = 264 h = 10.0 



a 


N 2 2 


M* 

z 2 


% deviation 


0.20 


0.5650 


0.5278 


6.58 


0.40 


0.7071 


0.6919 


2.15 


0.60 


0.6968 


0.6920 


0.69 


0.80 


0.6250 


0.6071 


2.86 


1.00 


0.5356 


0.5202 


2.87 


1.20 


0.4525 


0.4388 


3.02 


1.40 


0.3799 


0.3679 


3.15 


1.60 


0.3176 


0.3075 


3.18 


1.80 


0.2654 


0.2555 


3.73 


2.00 


0.2211 


0.2304 


-4.20 
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C. HASKIND'S RELATIONS AND ENERGY CHECK 



As previously indicated, important checks on the numerical 
technique and the final results can be made by application of 
Haskind's relations and the energy check which are developed 
in the Appendix. These checks were applied to all computer 
runs made and the percent difference, using the surface 
integration of the pressure results as the standard, is 
plotted in Figs. 19 and 20. 

In general, it was noted that dependence of the results 
on the grid size, N, is much greater for the forces obtained 
by a surface integration of the pressure when calculated by 
use of Haskind's relations. It can be seen from Figs. 19 
and 20 that as frequency decreases, the percent deviation 
decreases. At large values of a the error becomes quite 
large. The effect of decreasing the grid size or increasing 
the depth is also clearly shown. For the 120 point grid the 
errors are much larger than for the 240 point grid in the 
case of the h = 2.0 circular cylinder. Similarly the errors 
for the 120 point grid are much smaller for h = 6.0 than 
for h = 2.0. 
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Figure 19: PERCENT DEVIATION OF HASKIND'S RELATIONS FROM 

DIFFRACTION THEORY FOR A VERTICAL CIRCULAR CYLINDER 



PERCENT DEVIATION 




Figure 20: PERCENT DEVIATION OF HASKIND'S RELATIONS 

FROM DIFFRACTION THEORY FOR A SPHEROID 
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V. CONCLUSIONS 



A computer program based on linear wave theory has been 
developed to calculate wave forces acting on submerged 
objects of arbitrary shape. In order to reduce the computa- 
tion time, the restriction of symmetry with respect to two 
planes passing through the center of the object was 
incorporated . 

On the basis of the results presented in Section IV, 
the following conclusions appear justified: 

1. All of the checks and comparisons were successful 
and therefore it appears that the present method yeilds 
accurate results up to approximately a = 3.0 depending on 
the grid size. 

2. The horizontal force and moment coefficients in all 
cases increase with the wavelength parameter a, at first, 
reach a peak and then decrease as a increases. 

3. The vertical force coefficient for a floating 
spheroid starts with- a value of it, corresponding to zero 
size to wavelength ratio, and decreases rapidly at first 
and later slowly as a increases. 

4. Haskind’s relations and the energy check give 
valuable checks on the accuracy and show the effect of 
increasing the grid size or the depth. 
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APPENDIX 



HASKIND'S RELATIONS AND ENERGY CHECK 
A. HASKIND'S RELATIONS 

The dimensionless wave force and moment components, F^(t), 
on an object of arbitrary shape may be obtained by substituting 
for P' from (16) and (17) in (32) and simplifying. For the 
diffraction problem, this gives 



The scatter potential can be eliminated from (A.l) through 
the application of Green's reciprocal theorem, applied to G 
and u. as 



-G(x,y ,z;C,n,?)g^-(x,y,z) ]ds, 

in the region R, as shown in Fig. 21. Since an incompressible 
fluid was assumed, both u^ and u^ satisfy the Laplacian in 
in the region R, Eq . (18a), 




] , i = 1, 2, . . .6 (A.l) 



s 



where u^ denotes the radiation potentials. 
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Figure 21: REGION OF APPLICATION OF GREEN'S THEOREM 



(A. 2) 



II 



(u 



S+S f +S b +S ~ 



3u, 

i3n 



3u . 



-u 



73n 



i)ds = 0 



The contributions to the surface integral in (A. 2) from the 
surfaces S^, S u and S vanish since u 7 and u. satisfy the 
free surface and bottom boundary conditions and the radiation 
condition at R = °°. Also, from (18d) in conjunction with (14) 
and (15), (A. 2) yields 



rr 3u . rr 3u. 3u n 

/ / (u 7 + u ) 7T — — ds = / / (u - 5 — — -u. w——) ds . 
JJ v 7 o 7 3n JJ K o3n l 3n 7 



(A. 3) 



Again using Green's 
hand side of (A. 3) 

c r ■ 

JJ ( V u o 5 liT 

s 



theorem, the surface integral on the right 
may be replaced, thus giving 




(A. 4) 



By using the radiation condition (18e) and since x = r cos 
it can be shown that 





3u 

O', 



ds 




u u - [ iar (1-cos 
o 1 L v 

s 

oo 



e) 




ds . 



(A. 5) 



Substituting for u q and u^ from (14) and (15) and then 
simplifying 
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1 




a cosh 



1 



(ah) 



r sinh(2ah) h 1 iar 
L 4aE 2 J e 



2 TT 



J D.(0)r" 1/2 [l/2-iar(l - cos 9)] e iar cos 6 d0 



(A. 6) 



i ar cos 0 

Inasmuch as r 00 and e is a fluctuating quantity, 

the real part of (A. 6) equals zero. Thus 

r 2Tr 

- J D^(0)iar ^ 2 (1 - cos 0) e aar cos ® d0 



-1/2 iar 
= -iar ' e 



f 2 " 

J D i 



(0) (1 - cos 0) e 



-iar(l - cos 0) dQ (A. 7) 



Since ar >> 1, the solution to the integral on the right 
hand side of (A. 7) may be obtained by the method of stationary 
phase (Stoker [11]). Hence 




(0 ) (1 - cos 0) e 



- iar (1 



ar -> 00 



cos 0) d0 



(A. 8) 



= 2 ( 





D ± (tt ) 



and after simplifying 
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II 

c 



00 



Su . ; 3u 

(u o^r" u i ds 



2TT 1/2 
i (— ) 



,2 

cosh 



(ah) 



r sinh(2ah) il n , . 

[ ~ + h ] e 4 (tt) 



(A. 9) 



From (A.l), (A. 3), (A. 4) and (A. 9), simplification gives 



(A. 10) 



1/2 

F. (t) = Re{-i Ills’ [SinhUahi , h] e i( r »t)„ (l)) 

cosh z (ah) za 1 

and taking the modulus on both sides gives 
1/2 

C = 111 il f sinM2ah) + h ] | D (tt)| , i = i, 2, ...6 (A. 11) 

1 cosh Z (ah) Za 1 

In actuality, (A. 11) represents six different equations 

which are known as "Haskind's relations." They relate the 
t h 

i L component wave force (or moment) coefficient to the waves 

t h. 

produced at infinity by the object oscillating in the i mode. 



B. ENERGY CHECK 

t h 

Equation (31) represents the l component of the dynamic 

t h 

force (or moment) arising due to the j mode of oscillation 
of the object. In the case of the radiation problem, they 
may each by expressed as the sum of two components, one in 
phase with the acceleration and the other in phase with the 
velocity of the body in the form 



F. . (t) = - M..aX.(t) - N..aX.(t) , i , j = 1, 2, ...6 
ij v ' ij J ij 3 



(A. 12) 



where kb . represents the added mass (or added moment of 
inertia) and is defined, when compared to (34), (36) and (37), 
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as 



M. . 
ij 



-3 f C i 1 , 2 , 3 

= -pa Re [ J I Uj(x,y,z)h i (x,y,z)ds] , j = 1} 2, ...6 CA.13a) 



for a force and for a moment as 
M 



-4 f r i = 4, 5, 6 

ij ‘ - pa Re (/ J u j C x .y. z )h i (x,y,z)ds] , i . 1( 2 6> (A. 



13b) 



and denotes the damping coefficients (in relation to 

energy transport) and is defined, when compared to (34), (36) 

and (37) , as 



_3 r r i-1,2,3 

N i j = -paa Im[ / I (x,y,z)h i (x,y, z)ds] , . = 1 2 ,...6 ^ A - 14a ) 



for a force and 



N. 



-4 f f i = 4, 5, 6 

j = -paa Im[J I (x,y , z)h i (x,y , z)ds] , . = 1 , 2 ,...6 f A - 14b ) 



for a moment. 

til 

The dimensional force (or moment) component, in the i 
direction, on the oscillating body due to oscillations in the 
i*"* 1 mode may then be obtained from (A. 12) by equating j to i. 
Since the forces and moments exerted by the body on the fluid 
are equal and opposite to these, the energy transmitted by 
the object to the surrounding fluid during one period is given 
by 

T 

E. = -J (1 or a) F ii (t) L(t) dt (A. 15) 

O 

where a is used for a force (i = 1, 2, 3) and the factor 1 for 
a moment (i = 4, 5, 6). After carrying out the integration. 
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(A. 15) may be rewritten as 



E i = 7TpCT 2 a 5 X? 2 N ii} i = 1, 2, ...6 (A. 16) 

where denotes the dimensionless damping coefficient which 

is defined for a force as 

N.. i = 1, 2, 3 

N = ± - 1 . 

ij ’ j = 1, 2, . . .6 

J paa J ’ * 

and for a moment as 



N. . 
ij 



N. . 
ij 



paa 




5, 6 



2 , . . . 6 . 



If p? is the amplitude of the outgoing progressive waves 

til 

at S ro due to oscillation of the body in the l mode, then 
this amplitude can be related to the asymptotic velocity 
potential, u^, by using the dynamic free surface boundary 
condition and 



p?(r,0) = aX? tanh(kh) a |u^(r,G,0)| . (A. 17) 

The total energy transmitted over one period can now be 
written 



E. 

l 



Lim pL ri A 2kh 
4 ^ 

r-*- 00 sinh(2kh) 



a a x o 2 f 

~g i 



2tt 

| Uj (r , 0 ,0 ) | ^d0 , i = 1, 2, ...6. (A. 18) 

Now, by comparing (A. 16) and (A, 18), the damping coeffi- 
cients may be written in terms of the asymptotic velocity 
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potential , , .as 



N. . 
11 



L lm 1, i * ■» . r i 2a h -i 

^ t tanh(ah)[l + — UVo - ] 

r-*-°° 2 v J L sinhC2ah) J 



/ 



2 7T 



|u i (r,e,0) | 2 d0 , i = 1, 2, .. .6. 



(A. 19) 



C. NUMERICAL EVALUATION OF HASKIND'S RELATIONS AND ENERGY CHECK 
The right hand sides of (A. 11) and (A. 19), involving the 
asymptotic velocity potentials, are solved by the use of the 
distribution functions f., obtained from the numerical solu- 
tion of the radiation problem, in conjunction with the series 
form of Green's function (24). 

For r-^->°° and y = 0, considerable simplification is possible 
in the asymptotic form of G and may be written as 

~ r 2 2 . 

G (x , 0 , z ; £ , n , c) = - ? V -* cosh(ah) cosh [a(n + h) ] 

an- v n+v 



2 giCar^iT/4) 

O Y* * * 



v 7rar 

Thus the asymptotic potential, u^, is given by 



(A. 20) 



u. ( r , e. 0) . - hi-) 172 C v2 - a , 2 ) co . ; h ( a !l ) 

1 2 1,3 a 2 h-v 2 h+v 



If 



_ 1/2 ■ r IT. 

f^(£, n, O cosh [a(n + h)]r^ e 1 ^ 1 *! T-^ds . (A. 21) 



Since it is large, r^ may be replaced by r and moved outside 
the integral. Further, as shown in Fig. 22, for large r 

r^ = r-p^ cos(3-0) 
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Figure 22: PLAN VIEW OF THE OBJECT 



2 2 

where = (£ + £ ) , and ik may be rewritten as 



- 1/2 r ,2 2 



(v"-a")cosh(ah) _1 ^ 2 i(ar-T-) 

^ 277 4 

a h-v h+v 



(r , 0, 0 ) = -i (2 a) 

JJ £ i (?, n, c) cosh [a (ri+h) ] e la Pi cos (£" e ) ^s 4 

s 

Substituting (A. 22) into (A. 19) and simplifying gives the 
result 



(A. 22) 



, 2 2 
N . . = -i— [ 3 ' v 

1 1 /IT 1 



ii 2 tt l 2, 27 ~ , 

a h-v h+v o 



] / iff f 1 U, n, ? )cosh[a(n+h)] 



1/ 2s 



(e~ ia P! C0S C3 - 0 ) + e -iap 1 cos($+0)-j ds 



d6 . 



(A. 23) 



Once u^(r, 6, 0) is known, D^(tt) may be obtained by using 
the radiation condition (18e). Thus 

1/2 . 

D^(it) = u d (r,Tr,0)r e lar (A. 24) 

and after simplifying by means of trigonometric identifies, 

(A. 11) becomes 

C i ° c o~ sh (ah) JJ f i (C,v,Ocosh[a(n*h)]e' ia,r l cos(6 ' ,r) ds, (A. 25) 

l/2s 

where i = 1, 2, ♦•6 since for an axisymmetric body, the only 
cases of interest are those pertaining to i = 1 , 2, and 6. 

A comparison of the results obtained from (A. 23) and (A. 25) 
with those calculated by integration of the pressure over the 
submerged surface, Eqs . (37) and (38), are a valuable self- 

check on the accuracy of the numerical results. 
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